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I .   INTRODUCTION 

During  the  past  few  years,  there  have  been  many  advancements  in  the 
numerical  techniques  for  predicting  the  behavior  of  complex  fluid  flows. 
For  example,  several  computer  models  have  been  developed  by  Gosman, 
Spalding  and  others  [1,2,3]  which  use  the  mass,  momentum  and  energy 
conservation  equations  reduced  to  finite  difference,  nonlinear  algebraic 
form.   The  development  of  reliable  computer  programs  of  this  type 
greatly  benefits  engineering  analysis  in  such  widely  varying  fields  as 
meteorology,  aerodynamics  and  gasdynamics. 

The  earlier  two-dimensional  computer  codes  were  based  on  vorticity 
(oj)  and  stream  function  (4>)  [1,2,4].   This  form  of  the  governing  equa- 
tions eliminates  pressure  and  velocity  from  immediate  consideration. 
Pressure  is  normally  calculated  only  after  a  converged  solution  is 
obtained.   This  technique  has  several  inherent  disadvantages: 

1.  It  results  in  large  errors  in  the  predicted  pressure  distribu- 
tions in  all  but  quiescent  flow  regions  due  to  the  higher  order  depen- 
dence of  the  pressure  gradient  on  stream  function  [5]. 

2.  It  is  usually  restricted  to  constant  density  flows  or  to  flows 
in  which  density  varies  only  with  temperature  [3,5]. 

3.  The  boundary  conditions  are  difficult  to  specify  [3,5]. 

k.      Considerable  difficulty  has  been  experienced  in  arriving  at 
converged  solutions,  especially  for  nonuniformly  spaced  grids  and  high 
flow  rates  [2,4,5]. 

5.   The  ip-o)  model  is  not  easily  extended  to  three  dimensional 
flows  [3]. 


To  overcome  these  difficulties,  emphasis  has  been  placed  on  developing 
computer  codes  based  on  velocity  and  pressure,  the  primitive  variables. 

A  major  problem  with  any  new  computer  model  is  model  validation. 
The  difficulties  of  collecting  accurate  empirical  data  are  multiplied 
when  investigating  three  dimensional  and/or  reacting  flows.   In  addition, 
many  variables  within  these  flows  are  not  readily  measurable  (turbulence 
intensi  ties ,  etc.) . 

An  effort  to  utilize  elliptic  computer  models  which  can  handle 
turbulent,  reacting,  variable  density  flows  at  high  subsonic  and  sonic 
velocities  has  been  underway  at  the  Naval  Postgraduate  School  for  several 
years.   Two  specific  areas  which  have  been  investigated  are   flows  in  a 
turbojet  test  cell  and  in  the  combustion  environment  of  a  solid  fuel 
ramjet. 

A  solid  fuel  ramjet  (SFRJ)  most  often  consists  of  a  solid  fuel  grain 
which  provides  the  walls  for  the  combustion  chamber  [h] .      Located  at  the 
air  inlet  end  of  the  combustor  is  a  sudden  expansion  or  other  type  of 
flame  stabilization  device.   The  opposite  end,  downstream  of  the  fuel 
grain,  may  also  incorporate  a  sudden  expansion  aft  mixing  chamber.   The 
primary  combustion  region  contains  a  turbulent  diffusion  flame  which  emanates 
from  the  forward  recirculation  zone  and  remains  within  the  developing 
boundary  layer.   The  aft  mixing  region  may  incorporate  some  means  of 
injecting  air  (bypass  air)  in  order  to  complete  the  consumption  of  the 
fuel  which  exits  the  aft  end  of  the  fuel  grain.   Mixing  chamber  and  inlet 
design  variables,  fuel  grain  design  and  fuel  properties  make  a  wide  variety 
of  performance  characteristics  available. 


The  possibility  of  incorporating  this  type  of  propulsion  device  into  a 
future  medium  or  long-range  tactical  weapon  system  coupled  with  the  ex- 
pense of  testing  each  new  design,  makes  the  continued  development  of 
reliable  computer  models  highly  desirable.   The  model  could  be  used  to 
predict  the  effects  of  fuel  properties  and  to  inexpensively  evaluate 
different  geometries  and  operating  conditions.   In  addition,  a  three 
dimensional  code  would  allow  modeling  discrete  air  injection  into  the  aft 
mixing  region.   The  latter  technique  can  substantially  increase  combus- 
tion efficiency  and  allowable  fuel  loading. 

Previous  work  at  the  Naval  Postgraduate  School  has  been  directed 
toward  improvement  of  the  quantitative  accuracy  of  the  ip-oj  model  and 
toward  validation  of  that  model  [k] .      Reasonable  agreement  with  empirical 
data  has  been  obtained.   However,  as  previously  stated,  the  ip— 0)  model 
does  not  predict  accurate  pressure  distributions  and  numerical  difficul- 
ties prevented  modeling  the  aft  mixing  chamber. 

The  purpose  of  this  investigation  was  to  adapt  and  validate  a  primi- 
tive variable,  two-dimensional,  finite  difference  computer  code  which 
models  the  flow  within  a  solid  fuel  ramjet. 


I  I  .   MODEL  OVERVIEW 

A.  INTRODUCTION 

The  computer  model  used  in  this  study  was  adapted  from  the  CHAMPION 
2/E/FIX  computer  program  developed  by  Pun  and  Spalding  [6].   CHAMPION  is 
a  TWO-dimensional  Elliptic,  FIXed  grid  computer  program  which  provides 
a  solution  of  the  conservation  equations  for  recirculating  flows  in 
finite  difference  form. 

B.  ASSUMPTIONS 

The   flow  was   assumed   to  be   steady,    two-dimensional    and  subsonic. 
For  simplicity   the  value  of  specific  heat    (C   )    was   assumed   to  be  constant 
although    its   dependence  on   temperature   and/or  composition   could  easily 
be   included. 

A  modified   Jones-Launder    [6,7,8,9]    two  parameter   turbulence  model 
was    incorporated   to   calculate   the  effective  viscosity.      It   uses    five 
empirical    constants    (Table    I)    and    requires    that   two  additional    variables, 
turbulencekineti c  energy    (K)    and  turbulencedissipation    rate    (e) ,    be 
evaluated.      Effective   viscosity  was   calculated   using   the   formulas: 

^eff=Ulam+^t  (l) 

where 

ut  =   CD  p   K2/£  (2) 


cl 

C2 

CD 

ak,eff 

e  ,eff 

1.43 

1.92 

0.09 

1.0 

1.3 

TABLE  I.   K-e  TURBULENCE  MODEL  EMPIRICAL  CONSTANTS 


For  reacting  flows,  the  four  species,  oxygen,  nitrogen,  fuel,  and 
products,  were  considered.  Simple,  one-step,  infinitely  fast  kinetics 
were  assumed  in  which  a  fuel  combines  with  an  oxidant  to  form  a  single 
product  without  intermediaries  [^,10]. 

1  kg  fuel  +  i  kg  oxidizer     ■*■   (1  +  i)  kg  products 

Fuel  and  oxygen,  therefore,  could  not  exist  simul taneoulsy  and  the  com- 
bustion process  was  mixing  limited.   In  addition,  it  was  assumed  that  no 
oxygen  existed  at  the  fuel  surface  and  that  surface  was  isothermal.   The 
turbulent  Prandtl  and  Schmidt  numbers  were  taken  equal  to  unity  and, 
therefore,  the  turbulent  Lewis  number  was  unity.   The  laminar  Prandtl 
number  was  also  taken  to  be  unity. 

C.   GOVERNING  EQUATIONS 

The   conservation   equations    for  axi-symmetri cal    flows  with   no   tangential 
variations    can   be  put    into   the  general    form    [6]: 

f_(pu,)  .I§_<prv,)  -f^fi)  -  l-f^fi)    -  s+  (3) 

convection  terms        diffusion  terms  source  terms 

where  cj>  stands  for  the  dependent  variable  (u,v,k,£,h,  etc..)  being  con- 
sidered (<j>  =  1  for  the  continuity  equation),  T      is  the  appropriate  effec- 
tive exchange  coefficient  for  turbulent  flow  and  S,  is  the  "source  term" 
(Table  II).   The  energy  equation  in  terms  of  stagnation  enthalpy  has  no 
source  terms  since  the  turbulent  Prandtl  and  Schmidt  numbers  were  chosen 
as  unity  and  radiative  transport  was  neglected  [l,3l-   Thus  the 
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stagnation  enthalpy  is  given  by: 


h  =  h  +  (u2  +  v2)/2  +  K  (k) 


where  for  non- reacting  flows 


h  =  Cp  T  (5) 


and  for  reacting  flows: 


h  =   m   AH/i  +  C  (T  -  T  ,)  (6) 

ox        p      ref 


The  calculation  of  temperature  was  made  using  equations  (b) ,    (5) 
and  (6).   Density  was  calculated  from  the  perfect  gas  law: 

P  =  P/RT  (7) 


where, 


for  non-reacting   flows:    R  =  constant 
for   reacting    flows:      R  =   R/M 

1/M  =  m,   /Mr     T  mox     ox  N2      N2  pr     pr 

ru      ~u 


(8) 


(9) 


For  modeling  reacting  flows,  two  additional  quantities,  mN_  and 

v  =  mr  -  m  /i,  were  evaluated.   Each  of  these  properties  as  well  as 
A    fu    ox 

stagnation  enthalpy  have  identical  governing  differential  equations 
(equation  (3)  with  no  source  terms).   In  appropriate  dimens ion  less  form 
they  also  have  identical  boundary  conditions.   Thus,  only  one  of  the 
equations  had  to  be  solved.   The  dimens ionless  form  selected  for  each 
property  was : 


H  =  (h,n  -  h)/(h,n  "  hf)  (10) 


"W  =  (mN2.   "  m  N2)/(mN2.    ■  mN2,  )  (ll) 

in  in       rg 

x  -  (x  -  xin)/(xfg  -:xin)  (12) 


In  this  study,  stagnation  enthalpy  was  calculated.   H  was  then  formed 

using  equation  (10).  Since  H  =  m  2  =  x  at  all  points  in  the  flow  field, 

mN_  and  x  could  be  calculated  using  equations  (11)  and  (12).   The  mass 

fractions  of  fuel,  oxygen,  and  products  (mc    ,  m   ,  m  )  were  found  from 

73  r  fu*     ox'      pr 

the  equations: 

for  y  >  0;  mc     =  x»      m       =  0 

A  —  fu        A         ox 


(13) 


for  x  K  0"»  m£     -  °>     m       =  "Xi 

A  fu  ox  A 


m        =l-ni-m        ~mr      -   m.10  (14) 

pr  pr  ox  fu  N2 

D.      CONSERVATION   OF  MASS 

On  each    radial    line   the  mass    flow    rate  was    calculated   using    the    local 
density.      The  error    in   mass    flow    (compared    to   the  summation   of  "mass- in" 
at   all    upstream  boundaries)    was    used    to   uniformly   adjust    the  axial    veloci- 
ty over   the  entire    line.      This   process   ensured   that   overall    continuity  was 
satisfied  on    the    line.      The   pressure   at   all    downstream   locations  was    then 
adjusted   to  approximately   correct   for    the  momentum    imbalance   created   by 
the   uniform  axial    velocity   adjustment.      A  "pressure   correction"  equation 
was    then   solved   for  each   cell    on    the    line.      Local    cell    velocity    (axial 
and    radial)    and  pressure  were   then   adjusted   to  satisfy   cell-wise   conti- 
nuity.     The  details   of   this   procedure   are   presented    in    reference  6. 


E.      BOUNDARY   CONDITIONS 

1 .  Introduction 

Fixed  boundary  conditions  were  specified  at  the  desired  or  ex- 
perimentally determined  values.   Specified  gradient  boundary  conditions 
were  handled  by  setting  the  appropriate  convection/diffusion  coefficient 
to  zero  in  the  finite  difference  equation  ("breaking  the  link")  and  then 
entering  the  appropriate  gradient  through  linearized  "false"  source  terms 
[6].   The  geometry  and  appropriate  boundary  conditions  are  summarized  in 
figure  1 . 

2.  Inlet 

Although   not   a   computer  program   limitation,    "plug   flow"  was   assumed 
at   the   inlet.      Turbulence  kinet[c  energy  was   selected   to  be   uniform  with  a 
value  which   corresponded   to   the  approximate   turbulence    intensity  of  the 
inlet   flow. 

3.  Axis    of   Symmetry   and   Exit   Plane 

Radial    and   axial    gradients  were  set  equal    to   zero  on    the   center 
line  and  exit    respectively.      The    radial    flow   velocity  was  equated    to     zero. 
k.      Sol  id   Boundaries 

All    non-reacting   solid  boundaries  were   considered  adiabatic  with 
both   velocity   components   equal    to  zero    ("no  slip"   condition). 

For  simplicity,    a   two  part   boundary    layer  was    used.      The  border 
between   the    laminar  sublayer  and   the   turbulent    layer  was    taken   at 

y       =11.5    [6].      y       was   evaluated   at  each   near  wall    node    (p) , 

<-(<^,am>(vp»'/2  <15» 


10 


where,    for  y        >    11.5 

Tw=CD,/2PKp  ("«> 

xw  was   assumed  uniform  from  the  wall    to  the  near  wall  grid  point.      Thus, 

Yp+  "  'd'''  P  Kp'/2  «/u  ,„  C7) 


If  y       >_  11.5,    the  wall    shear  stress    (x  )   was   calculated  using    the   formula 


Tw       CD         p   Kp  "  p  CD         Kp        (u/u  } 

,  r    1A  n  „    1/2,.    ,_    ..    ]/k      1/2,  v 

=  K  CD  p   up   Kp        /^^P^d       Kp        ^1am) 


(18) 


where 


u+-lln(Eyn+)  (19) 

Wall    shear   stress  was   evaluated   for  y        <    11.5    from   the  formula: 

P 


Tw  -  ^lam  Up/5  <20> 


Due   to   the  steep  gradients   of  properties    in    turbulent    flows   near 
solid  boundaries,    the  source   terms    for  K  and  £   at   near  wall    nodes    were 
expressed    in    terms   of   the  wall    shear  stress    [1,6].      T     also   provides    the 
boundary   condition    for   the  u  and   v  equations.      In    the   following  equation 
for   turbulence  dissipation    rate    (e)    at  a  near  wall    node    (p) ,    the   length 
scale    is   presumed  proportional    to   the  distance   from  the  wall    (6). 

£     =  Cn3A  K  3/2/k6   =  K  3/2/2.436  (21) 

p  D  p  p 
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It  was  found,  as  was  previously  found  by  Netzer  [k] ,    that  (when  using  the 

sudden  expansion  geometry  in  reacting  flows)  the  near  wall  dissipation 

3/2 
had  to  be  increased  on  the  step  face  (e  =  K    /0.46)  and  that  the  grid 

P    P 

spacing  adjacent  to  the  fuel  surface  had  to  be  fine  (y   <  11.5)  in 
order  to  obtain  a  temperature  distribution  in  qualitative  agreement  with 
experiment.   Equation  (20)  implies  that  the  wall  shear  stress  is  calcula- 
ted assuming  a  linear  velocity  profile  when  y   <  11.5.   A  near-wall 
grid  point,  therefore,  can  lie  within  the  laminar  sublayer,  but  the 

source  terms  for  K  and  z    imply  that  u  ££/ii,    is  much  greater  than  one 

r    '  eff      lam  3 

[7,8].      This    fact  precludes   y        from  being   significantly    less    than    11.5. 

For   reacting    flows,    the  boundary   conditions    for   the  dimensionless 
properties    (equations    (10),    (11)    and    (12))    were   zero  at   the   inlet   and 
unity  "deep"    in    the   fuel    grain    (fg).      These  properties   were   considered 
to  have   zero   gradients   on   non-reacting   surfaces. 

The  assumptions   employed   for   reacting   flows    (unity   turbulent   Prandtl 
and   Schmidt   numbers,    simple   chemical    reaction,    constant   specific  heat   and 
stagnation  enthalpy   defined    in   equations    (k)    and    (6))    result    in   a  general 
boundary   condition   for  all    "conserved"   properties    (<J)  )[10]   on   a  surface 
which   has  mass    transfer, 

<       =    (I\  -K=r)u   /(<L        -   4>        )  (22) 

bw  <f>   3r      bw     ^        ^cfg 

where  <$>      represents   h,    m..0   or  y   =  mr     -  m     /i  . 
c       r  '      N2  A  fu  ox 

A  mass    transfer  cnnductance    (g)    is   often   defined   such    that, 

3<f> 

(i\  *-£)u     »  g(<f>    "  <j»      s  (23) 

d>   9r     bw        aNTc  c,    ) 

T  oo  bw 


where  $        is  defined  as  the  free  stream  value.   For  this  application,  0 

Coo 

was    taken    to  be   the    local    near  wall    value  <|> 


OO 


CP 


Substituting  equation  (23)  into  equation  (22)  yields 
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p     bw     bw     fg 
=   9  BP  (25) 

where  BP  represents  the  mass  transfer  (or  "blowing")  parameter. 

Without  mass  transfer  the  wall  heat  flux  (q  ")  can  be  defined  in 

^w 

terms  of  the  conditions  at  the  near  wall  node. 


4« "  ■  Hw  -  V  -  -  <M?)«  <26> 

p     y  p 

where  h   is  the  enthalpy  of  the  wall  and  h   is  the  heat  transfer  conduc- 
w  r/ 

tance. 

With  9  -  h  in  equation  (23), 

g  -  (£-!£)  /(hD  -  h  J  (27) 

C   or  w   P    w 

Substituting  equation  (26)  into  equation  (27)  yields: 

g  =  h/Z  (28) 


or 


thus, 


g/(pu)p  =  fe/[(pu)p  Cp]   E  St  (29) 


g  =  (pu)   St  (30) 


From  Reynolds  Analogy  with  unity  Prandtl  number, 


St  =  C,/2  =  t  /(pu2)  (3D 

r      w      p 


where  C^  is  the  local  friction  coefficient. 
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Combining  equations  (30)  and  (31)  yields 


9  =  t  /u  (32) 

3    w  p  v 


Using   the  Couette   flow  approximation   for  the  boundary    layer  behavior 
with   mass    transfer    [10], 

g  =   g*    ln(l    +   BP)/BP  (33) 

where 

9*  =   I'mgp^tg)  (3*) 

In  this  application,  BP  was  evaluated  from  the  solution  of  the 
energy  equation  using, 

BP  -  (hp  "  \J'\„  -   hf„)  (35) 

The  wall    shear   stress  was    calculated  using   equation    (18)    or  equation 
(20)    and  modified  with   equation    (33). 

Tbw=   Tw    ln(1    +   BP)/BP  (36) 

where  t   is  the  wall  shear  stress  without  wall  mass  addition, 
w 

The  mass  transfer  conductance  (g)  was  found  using  equation  (32).   The 

wall  mass  flux  was  then  evaluated  using  equation  (25). 

The  wall  heat  flux  (q  ")  on  all  solid  isothermal  boundaries  was  evalu- 

w 

ated  using  the  Reynolds  analogy: 

-q  "/(h   -  h  )  =  T  /u  (37) 

^w    p    w'    w  p 


}k 


Since   the  blowing    rates  were  small    for   the  solid   fuel    ramjet    (typi- 
cally,   BP  <   2.0),    K  and  e  were  evaluated   using   equation    (3)    and   the   terms 
presented    in   Table    II    which    incorporate   the  empirical    constants  of  Table 
I. 

Blowing   velocity    (v,     )    and   fuel    regression    rate   (RR)    were   calcula- 
ted  using   the   formulas: 

vbw  "  -V'^bw  <38) 

RR    =     mbw"/pfg  (39) 

F.      SOLUTION   PROCEDURE 

Five   variables    (u,v,K,£   and   H  or  h)    were   solved   using  equation    (3)    in 
finite  difference   form.      The   line  by    line   iterative   procedure  employed 
upwind   differencing   and   under   relaxation    to  promote   convergence    [6]. 
Pressure    (relative   to   a  selectable   position   and  magnitude  within    the 
grid)    was   obtained   from  the  mass    conservation    imposed  on   each    radial 
grid    line   and  on  each   nodal    control    volume   as   discussed  above.      Effective 
viscosity,    temperature   and   density  were   also  obtained   as    described   above. 
A  more  detailed  explanation  of   this   procedure   can   be   found   in    reference6. 
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III.   DISCUSSION  OF  RESULTS 

A.   INTRODUCTION 

The   purpose   of   this   study  was    to  develop   a   primitive   variable,    finite- 
difference   computer   program   that   could  be   used   to  determine   the   flow  within 
a   solid   fuel    ramjet   combustor  with   emphasis   on   the  aft  mixing   chamber.      The 
effects   of    inlet  mass    flow    rate   and    inlet   dump   area    ratio  on    the   flow   field 
were  examined.      As    previously   explained,    an   aft  mixing    region   allows    further 
combustion   aft  of   the   fuel    grain.      This   process   normally    increases    combustion 
efficiency.      Lowering    the    inlet   flow   rate    increases    the   fuel-air   ratio  within 
the   fuel    port.      Bypass   air  can    then    be    injected    into   the   aft   mixing    region. 
The    latter  procedure   can   be   used   to  appreciably    increase    fuel    loading.      Pre- 
vious work   at    the  Naval    Postgraduate   School     [4,11,12]    modeled   a   SFRJ  with   a 
computer  program  utilizing  i^-co  as   primary   variables.      Numerical    instabilities, 
however,    prevented   the   use  of   the  ifj-to  model    to   predict   the   flow    in    the   aft 
mixing    region.      The    results  of   that    investigation    and   some  empirical    data  were 
available   for   comparison  with    the  predictions    from   the  primitive  variable 
mode  1 . 

Several    factors  were   anticipated  which   could   contribute   to  differences 
between    the   predictions   of    the   two  models   and   the  empirical    data: 

a.  Some   of    the  experimental    data  were   measured    in    cold,    nonreacting 
f  1  ows . 

b.  The    incorporation  of  the  aft  mixing   chamber    into   the  primitive 
variable  model    could    influeoce   the   flow   upstream   in    the   combustion   chamber. 

c.  In    the  ip-co  model,    a  wall    value  of   turbulence  kinetic  energy    (K) 

was   specified   through   a   slip   factor  such    that   K     =    (-1.0   or   -0.39)*K   . 

w  p 
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depending   on    the  magnitude   of   the   turbulent    Reynolds    Number.      In    the 
u-v-p  model    the  boundary   condition    for   K  at   the  near  wall    node    (p)   was 
specified    in   terms   of   the  wall    shear  stress.    In   addition,    in    the   primi- 
tive variable  model,    the  boundary   condition   for  stagnation   enthalpy   at 
the  near  wall    node  was   made  a   function  of  wall    shear  stress    through    the 
Reynolds  Analogy.      These   factors   affect  heat   flux   to   the  wal 1 ,    and, 

therefore,    the   fuel    regression    rate. 

d.      The   u-v-p  model    incorporated   a  23  by  21    grid    in    the   fuel 
port  while   the  iJ/-to  model    utilized   a    17  by  25   grid.      In    reality   the 

heat  of  vaporization  of   the   fuel    is   a   fixed  quantity   and,    if   converged 
solutions    are  obtained,    the  wall    heat    flux  should  not   depend   upon    the 
grid  spacing.      However,    it  has   been    found    [k]    that   the  heat   flux  to   the 
wall    (which    is    calculated   using   the  near-wall    grid  point)    is    a 
function   of   the  grid   distance   from  the  wall.      This    results    from   the 
assumed  behavior  of   the  variables   near  the  wall.      The   procedure  employed 
in    this   study  was    to  adjust   the  heat  of  vaporization    to  match    the  empiri- 
cal   fuel    regression    rate  at  one  air   flow   rate  and   to  use   that   value  for 
all    other   flow   rates.      If   the  model    is    realistic,    fuel    regression    rate 
should   then   vary  with   air   flow   rate    in   agreement  with   experiment. 
B.      REGRESSION    RATE 

Figure  2   shows    that   the   fuel    regression    rate   predictions   of  the   u-v-p 
and  iJj-u)  models   are  quite   similar.      Both   predict   the  peak    regression    rate 
upstream  of  experiment   and  have   similar   slopes.      This   early   peak    in    the 
regression    rate    results    from   the  model    predicting   a  shorter    reattachment 
length    than  was    found  experimentally    [4].      The   primitive   variable  model 
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predicted  higher   regression    rates   downstream  of   flow    reattachment    in 

better  agreement  with   experiment. 

Figure   3   shows    the  effects   of   increasing    inlet   air  mass    flux    (G  = 

m   .    /A,,     a  V.    )    on    fuel    regression    rate    (rf  ).      The   regression   profile 
a i r     rp  in  r u 

remained   the  same  and,    as   expected,    decreased  with   decreasing   G.      It 

has   been    found  experimentally      that   the   regression    rate  of  plexiglass 

varies   as    the  air  mass   flux   raised   to  a   constant   power    (r.     °°  G  )  .      Boaz 

fu 

and   Netzer    [11]    found    that    this    constant  was   equal    to   0.41   while  Mady , 
et   al    [12]    found    it   to  be  approximately  0.38.      For   the   three   test   cases 
of   this   study,    the   u-v-p   model    predicted   the   constant,    n,    to  be  between 
0.31    and  0.3^-      Thus,    the  primitive   variable  model    appears    to  correctly 
predict    the  nature  of  the   change    in   convective   heat    flux   to   the   fuel 
surface  with   air  flow    rate. 
C.      TURBULENCE    INTENSITY 

Figure   k   comparies    the   predicted   centerline    turbulence    intensity 
(assuming    isotropic   turbulence)    and  experimental    data   for  non-reacting 
flow.     The   primitive  variable   computer  model    slightly   underpredi cted   the 
peak   turbulence    intensity  while   the  ifj-co  model    overpredi  cted    it.      Both 
models   predicted   the   peak   occurring   downstream  of  experiment   and   both 
distributions    appear   to  approach   an    identical    asymptote   downstream. 
The   decrease    in    turbulence    intensity   predicted  by   the  ip-U)  model    near   the 
inlet    resulted   from  the  model   over-predicting   the  velocity    increase   as 
the   air  entered    the   combustor    [k] .      The   u-v-p  model   overcame   this    diffi- 
culty.     The   differences    in    the    results    from   the   two  computer  models   may 
result   from   the   differences    in    the  boundary   conditions   on    turbulence 
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kinetic  energy  in  the  combustor  and/or  to  the  effects  of  the  addition 
of  the  aft  mixing  chamber  on  the  upstream  flow.   It  should  also  be  noted 
that  the  experimental  data  used  in  this  comparison  were  obtained  in  a 
non- reacting  flow. 

Figure  5  shows  the  effect  of  decreasing  inlet  air  mass  flux  on 
turbulence  intensity.   As  anticipated,  the  peak  turbulence  intensity 
decreased  as  inlet  axial  velocity  decreased.   Each  test  condition,  however, 
converged  on  the  same  value  downstream.   Much  additional  experimental  work 
is  required  to  obtain  the  turbulence  intensities  in  reacting  flows;  only 
then  can  the  adequacy  of  the  K-e  turbulence  model  be  fully  evaluated. 

D.  PRESSURE 

Figure  6   shows    the  effect   of    inlet   velocity  on    the  axial    pressure 
distributions    for   the    three  primitive   variable    test  conditions    (Table    III). 
(The    radial    scale   has   been   expanded   to    illustrate   the   pressure   variations. 
The  maximum  pressure   variation    is    approximately    1.2   psi.)      The    radial 
location  of  these  distributions   is    given    as    a    fraction    of    the'  fuel   port 
radius    (Rf  ).      As   expected,    pressure    initially    increased   due   to  jet 
spreading.      This  was    followed   by  a   slight  pressure   drop   as    the   flow   ac- 
celerated  due   to  heat   addition   and  wall    friction.      The   final    pressure 
rise  was   due   to  jet   spreading    In    the   aft  mixing    region. 

E.  TEMPERATURE 

Figure   7   displays    radial    temperature   variations    in    the   combustor 
near  the  end  of   the   full    grain   and   at   about    1.5   aft  mixing    region    dia- 
meters   down    the  aft   chamber.      As    discussed   above,    fuel    flow    rate 
decreases    less    than    inlet   air   flow   rate    (r,     °°  G    ,    n   <    1 )  .      Therefore, 
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as   air   flow   rate    is    decreased,    the  overall    mixture    ratio   becomes   more 
fuel    rich,    and   the   developing   boundary    layer  and   the   fuel    layer  between 
the   diffusion    flame   and   the  wall    thicken    .      Thus,    as    shown    in   figure  7,    as 
the   inlet   velocity    (and,    therefore,    the   inlet   air  mass    flux)   was    de- 
creased,   the  maximum   temperature    (or  "flame")    in   the  combustor  moved   away 
from   the   fuel    grain   and   the  centerline  temperature    increased.      The 
maximum   temperature    in    the   aft  mixing   chamber  was    also   predicted   to  occur 
farther  from  the   top  wal 1 . 

Figure   8  shows   similar  data  predicted   by    the  ty-u)  model    slightly 
farther   upstream.      A  significant   difference   between    the   predictions   of 
the   two   computer  models  was    that    the  ip-oo  model    predicted   a  stronger 
dependence  of   the  peak   temperature    radial    location  on    the   inlet   air   ve- 
locity.     An   aft  mixing    region  was   not    incorporated    into   the  ip-oo  model. 
Therefore,    the   boundary    layer  continued    to  grow   and   the   point  of  maximum 
temperature   continued   to    recede   from   the   fuel    surface  with    increasing 
axial    distance   from  the    initial    reattachment   point.      The   aft  mixing    region 
of   the  u-v-p  model    caused   the  boundary    layer   thickness    (and,    therefore, 
the    location   of   the  peak   temperature)    to   become  approximately   constant    in 
the    latter  portion  of   the   combustion   chamber.      This  was    the  apparent 
cause  of   the  weaker  dependence  predicted   by    the  u-v-p  model    of  peak 
temperature    location   and  boundary    layer   thickness   on    inlet   air  mass    flux. 

Figure  9    is    an    illustration   of   the  predicted   combustion   behavior    in 
the  aft  mixing    region.      (The    radial    dimension   has   been    expanded   for   clarity.) 
Lines   of  maximum  temperature    (i.e.,    the   flame  sheet    location)    are   pre- 
sented as    a   function  of  fuel    grain    inlet   air  velocity.      It   should  be  noted 
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that   the  aft    recirculation  zone,  which    is   also  depicted  on   this   figure, 
was   predicted   to  be   fuel    rich   and   did   not  vary   appreciably    in   size  with 
changing    inlet  air  mass    flux.     As   discussed   above,    the   fuel    regression 
rate  decreased  more  slowly   than    inlet  air   flow   rate.      Thus,    as   air 
flux  was   decreased   the  mixture  entering   the   aft  chamber  became  more   fuel 
rich  and   the  thickness  of   the  fuel    layer  at   the  end  of  the  fuel    grain 
increased  slightly.      With   high   air  mass    flux   through    the   fuel    port   the 
mixture    ratio    is    fuel    lean.      The   flame   therefore   propagates    to   the  outer 
wall    of   the   aft  mixing   chamber.      This   condition   could   be  expected   to 
produce   a  high   combustion   efficiency.      ft   should  be  noted   that   an   ade- 
quate   length-to-diameter   ratio    is    required   to  allow   the   flame   to   spread 
to   the  wall.      This    ratio    is   apparently   a   function  of  the   fuel    port   to   aft 
mixing   chamber  cross-sectional    area   ratio-       As    the   air   flow   rate  was 
decreased   the  mixture    ratio  became   fuel    rich   and   the  flame   did  not    reach 
the  wall-      This  would   result    in   unburned   fuel    entering   the  nozzle   and  a 
lower  combustion  efficiency.      Figure    10   shows    the   predicted  effect    of   the 
fuel    inlet   dump  step   size  on   the   flame  behavior    in  the   aft  mixing   cham- 
ber.       The    recirculation   zone   did  not  change    in   size  since   the   air  mass 
fluxes  were    identical.      The   smaller   inlet   step   produced   a  slightly   higher 
fuel    regression    rate   and   therefore    required   a    longer  aft   mixing    region. 
These  predictions   might  be   used   as    a   first  approximation   for  predicting 
the  "best"  placement  of  bypass    air  dumps    in   the  aft  mixing    region.      To 
predict   an  optimum   location,    however,    the  primitive   variable  model   would 
have   to   be  expanded   to   three  dimensions. 
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F.   COMPUTER  RELATED  PROBLEMS 

As   has   been   discussed   previously,    in   order   to  obtain    results    that  were 
in   agreement  with  experiment,    the  grid   spacing   near   the   fuel    surface  was 
required   to  be   fine   and   the    length   scale  of   turbulence  was    decreased  on 
the   combustor  step    face.      Because   convergence  was   sensitive   to  the    length- 
to-width    ratio'of   individual    control    volumes,    the  small    radial    grid   spacing 
near   the   fuel    surface   forced   similar   fine   spacing    in    the   axial    direction 
downstream    in   the   aft  mixing    region.      A    length    to  width    ratio  of    less    than    ten 
to  one  was    required.      These   criteria   forced   the   use  of  a    large  number  of  cells, 
which    in    turn    required  a    large   amount  of   CPU   time.      A   typical    primitive 
variable   *+0   by   33   grid    required   75    to   80  minutes,  of   CPU   time  on   an    IBM   360-67 
computer   to   converge.      A   typical    iJj-co  model    with    a    17   by   25   grid    required    35 
to   kO   minutes   of   CPU    time.       It  must   be    remembered,    however,    that   numerical 
instabilities    prohibited   the  modeling   of   an   aft  mixing   chamber  with    the  i^-oo 
model . 

The   primitive   variable  model    demonstrated   some   convergence   difficulty    in 
the  aft    recirculation    region.      This   problem  seemed    to  be   associated  with 
the   continually    changing   velocity   profile   just   prior   to   the   aft  expansion 
(the  "inlet"    conditions    for   the   aft   mixing    chamber).      This   effect  was 
suppressed   by    sweeping   through   the  entire   flow   field   several 
times  with   only   a    few   traverses   on   each    line  and   then    increasing    the  number 
of    traverses   on    the    radial    grid    lines    in    the  aft  mixing    region   once    the 
combustor    flow  field   had  essentially    converged. 
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IV.   CONCLUSIONS 

In  general,  the  predicted  flow  fields  for  the  two  computer  models 
were  quite  similar  within  the  fuel  grain.   However,  the  presence  of  the 
aft  mixing  region  coupled  with  the  few  boundary  condition   differences 
previously  mentioned,  had  some  effect  on  the  flow  field  predictions.   The 
most  noticeable  of  these  was  the  decrease  in  dependence  of  the  boundary 
layer  thickness  and  the  maximum  temperature  radial  location  on  axial  in- 
let velocity.    As~antT  cipated,  the  primitive  variable  model  allowed  the 
prediction  of  the  flow  field  within  the  aft  mixing  region.   This  was  not 
possible  with  the  ijj-oc  model.   Many  additional  empirical  data  are  needed 
to  completely  assess  the  validity  of  the  primitive  variable  model  in 
predicting  the  flow  in  a  SFRJ. 
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